import numpy as np
import matplotlib.pyplot as plt
from math import *

s_D = 0.1e-2

lambd = 632.8*10**-9        # Láser He-Ne

# Rendija
D1 = 194.8e-2  #cm
d_min_1 = [0.85, 0.8, 0.8]  #cm
n_1 = [-2, -1, 1]

D2 = 342.0e-2  #m
d_min_2 = [1.35, 1.45, 1.4] #cm
n_2 = [-2, -1, 1]
#\Deltax=x_{n+1}-x_n=lambda*D/a

# Agujero
D_agujero_1 = 154.4   #cm
d_min_agujero_vertical_1 = [0.3, 0.3, 0.35] #cm
d_min_agujero_horizontal_1 = [0.6, 0.6, 0.6] #cm
n_agujero_1 = [-3, -2, 2]

D_agujero_2 = 342.8   #cm
d_min_agujero_vertical_2 = [0.7, 0.7, 0.75] #cm
d_min_agujero_horizontal_2 = [1.3, 1.3, 1.3] #cm
n_agujero_2 = [-2, -1, 1]

# Medidas directas cuadrado
tamano_cuadrado = [[210, 210, 220], [220, 230, 220]]    #[[medidas anchura], [medidas altura]]   en micrometros

# Medidas directas rendija
tamano_rendija = [170, 175, 170]    #[anchura]   en micrometros

# -----------------------------------------------------------------------------

# Rendija
Dx1, Dx2 = np.array(d_min_1)*1e-2, np.array(d_min_2)*1e-2   #Pasamos a metros

A1, A2 = np.mean(Dx1), np.mean(Dx2)
sigma_A1, sigma_A2 = np.std(Dx1, ddof=1), np.std(Dx2, ddof=1)

a1, a2 = lambd * D1 / A1, lambd * D2 / A2
s_a1, s_a2 = sqrt(((lambd)/(A1))**2 *s_D**2 + ((lambd*D1)/(A1**2))**2 *sigma_A1**2), sqrt(((lambd)/(A2))**2 *s_D**2 + ((lambd*D2)/(A2**2))**2 *sigma_A2**2)

a, s_a = np.mean([a1, a2]), np.mean([s_a1, s_a2])

print("///Datos Rendija///" + "="*56)
print(f"A1 = {A1:.2e} +- {sigma_A1:.1e} | A2 = {A2:.2e} +- {sigma_A2:.1e}")
print(f"a1 = {a1:.2e} +- {s_a1:.1e} m | a2 = {a2:.2e} +- {s_a2:.1e} m\n\t\t\t\t\t\t-----> a = {a:.2e} +- {s_a:.1e} m")
print("="*75)

#Reconstruimos gráficas acumulando distancias. Es aproximado pero válido en la práctica
#Gráfica 1 (Set de datos 1, para D1)
dx1 = np.array(d_min_1)*1e-2

x1 = np.array([
    - (dx1[0] + dx1[1]),  # n = -2
    - dx1[1],             # n = -1
    + dx1[2]              # n = 1
])

n1 = np.array(n_1)

p, cov = np.polyfit(n1, x1, 1, cov=True)
A_fit, B_fit = p
sigma_A_fit = np.sqrt(cov[0,0])
sigma_B_fit = np.sqrt(cov[1,1])


a_fit1 = lambd * D1 / A_fit
sigma_a_fit1 = a_fit1 * sigma_A_fit / A_fit

print("\n///Ajuste lineal (D1)///" + "="*50)
print(f"A1 = {A_fit:.2e} ± {sigma_A_fit:.1e} m")
print(f"B1 = {B_fit:.2e} ± {sigma_B_fit:.1e} m")
print(f"a (desde ajuste) = {a_fit1:.2e} ± {sigma_a_fit1:.1e} m")
print("="*75)


n_fit = np.linspace(min(n1)-1, max(n1)+1, 100)

plt.figure()
plt.scatter(n1, x1, label="Datos experimentales")
plt.plot(n_fit, A_fit*n_fit + B_fit, label="Ajuste lineal")

plt.xlabel("n")
plt.ylabel("x (m)")
plt.title(f"Difracción por rendija (D = 194.8 cm)")
plt.legend()
plt.grid()

plt.show()

#Gráfica 2 (Set de datos 2, para D2)
dx2 = np.array(d_min_2)*1e-2

x2 = np.array([
    - (dx2[0] + dx2[1]),  # n = -2
    - dx2[1],             # n = -1
    + dx2[2]              # n = 1
])

n2 = np.array(n_2)

p, cov = np.polyfit(n2, x2, 1, cov=True)
A_fit, B_fit = p
sigma_A_fit = np.sqrt(cov[0,0])
sigma_B_fit = np.sqrt(cov[1,1])


a_fit2 = lambd * D2 / A_fit
sigma_a_fit2 = a_fit2 * sigma_A_fit / A_fit

print("\n///Ajuste lineal (D2)///" + "="*50)
print(f"A2 = {A_fit:.2e} ± {sigma_A_fit:.1e} m")
print(f"B2 = {B_fit:.2e} ± {sigma_B_fit:.1e} m")
print(f"a (desde ajuste) = {a_fit2:.2e} ± {sigma_a_fit2:.1e} m")
print("="*75)


n_fit = np.linspace(min(n2)-1, max(n2)+1, 100)

plt.figure()
plt.scatter(n2, x2, label="Datos experimentales")
plt.plot(n_fit, A_fit*n_fit + B_fit, label="Ajuste lineal")

plt.xlabel("n")
plt.ylabel("x (m)")
plt.title(f"Difracción por rendija (D = 343.0 cm)")
plt.legend()
plt.grid()

plt.show()

#Resultado de a según ajustes
a_fit, s_a_fit = np.mean([a_fit1, a_fit2]), np.mean([sigma_a_fit1, sigma_a_fit2])
print(f"\n///\t\t\t\t ----> a (desde ajuste) = {a_fit:.2e} +- {s_a_fit:.1e} m\n")
print("-"*75)
a_t = np.mean([a_fit1, a_fit2, a])
s_a_t = np.mean([sigma_a_fit1, sigma_a_fit2, s_a])
print(f"\n///\t\t\t\t ----> a (total) = {a_t:.2e} +- {s_a_t:.1e} m\n")
print("="*75)

# Medida de anchura de la rendija con el microscopio:    

a_media_microscopio = np.mean(tamano_rendija)*1e-6
s_microscopio = 1e-6
print("\n///Ancho Rendija por Microscopio///" + "="*40)
print(f"\n///\t\t\t\t ----> a (microscopio) = {a_media_microscopio:.2e} +- {s_microscopio:.1e} m\n")


# -----------------------------------------------------------------------------  
# =============================================================================
# RECTÁNGULO AGUJERO (Abertura rectangular)
print("="*75)
print("="*75)
print("\n///RECTÁNGULO///" + "="*59)

print("\nD1 = 154.4 cm ------------------------------------------------------------")

# Pasamos D a metros D1
D_rect = D_agujero_1 * 1e-2

# Pasamos a metros
dx_h = np.array(d_min_agujero_horizontal_1) * 1e-2
dx_v = np.array(d_min_agujero_vertical_1) * 1e-2

# -----------------------
# Cálculo con Δx D1

A_h = np.mean(dx_h)
A_v = np.mean(dx_v)

sigma_A_h = np.std(dx_h, ddof=1)
sigma_A_v = np.std(dx_v, ddof=1)

# Dimensiones
a_rect_1 = lambd * D_rect / A_h
b_rect_1 = lambd * D_rect / A_v

sigma_a_rect_1 = sqrt(((lambd)/(A_h))**2 * s_D**2 + ((lambd*D_rect)/(A_h**2))**2 * sigma_A_h**2)
sigma_b_rect_1 = sqrt(((lambd)/(A_v))**2 * s_D**2 + ((lambd*D_rect)/(A_v**2))**2 * sigma_A_v**2)

print("\n--- Método Δx ---")
print(f"A_h (para a) = {A_h:.2e} ± {sigma_A_h:.1e} m")
print(f"A_v (para b) = {A_v:.2e} ± {sigma_A_v:.1e} m")
print(f"a = {a_rect_1:.2e} ± {sigma_a_rect_1:.1e} m")
print(f"b = {b_rect_1:.2e} ± {sigma_b_rect_1:.1e} m")

# -----------------------
# RECONSTRUCCIÓN + AJUSTE (horizontal) D1

x_h = np.array([
    - (dx_h[0] + dx_h[1]),
    - dx_h[1],
    + dx_h[2]
])

n_h = np.array(n_agujero_1)

p_h, cov_h = np.polyfit(n_h, x_h, 1, cov=True)
A_h_fit, B_h_fit = p_h
sigma_A_h_fit = np.sqrt(cov_h[0,0])
sigma_B_h_fit = np.sqrt(cov_h[1,1])

a_rect_fit_1 = lambd * D_rect / A_h_fit
sigma_a_rect_fit_1 = a_rect_fit_1 * sigma_A_h_fit / A_h_fit

print("\n--- Ajuste horizontal ---")
print(f"A = {A_h_fit:.2e} ± {sigma_A_h_fit:.1e} m")
print(f"B = {B_h_fit:.2e} ± {sigma_B_h_fit:.1e} m")
print(f"a = {a_rect_fit_1:.2e} ± {sigma_a_rect_fit_1:.1e} m")

# Gráfica horizontal
n_fit = np.linspace(min(n_h)-1, max(n_h)+1, 100)

plt.figure()
plt.scatter(n_h, x_h, label="Datos experimentales")
plt.plot(n_fit, A_h_fit*n_fit + B_h_fit, label="Ajuste lineal")

plt.xlabel("n")
plt.ylabel("x (m)")
plt.title("Difracción por abertura rectangular (D = 154.4 cm) dirección horizontal")
plt.legend()
plt.grid()
plt.show()

# -----------------------
# RECONSTRUCCIÓN + AJUSTE (vertical) D1

x_v = np.array([
    - (dx_v[0] + dx_v[1]),
    - dx_v[1],
    + dx_v[2]
])

n_v = np.array(n_agujero_1)

p_v, cov_v = np.polyfit(n_v, x_v, 1, cov=True)
A_v_fit, B_v_fit = p_v
sigma_A_v_fit = np.sqrt(cov_v[0,0])
sigma_B_v_fit = np.sqrt(cov_v[1,1])

b_rect_fit_1 = lambd * D_rect / A_v_fit
sigma_b_rect_fit_1 = b_rect_fit_1 * sigma_A_v_fit / A_v_fit

print("\n--- Ajuste vertical ---")
print(f"A = {A_v_fit:.2e} ± {sigma_A_v_fit:.1e} m")
print(f"B = {B_v_fit:.2e} ± {sigma_B_v_fit:.1e} m")
print(f"b = {b_rect_fit_1:.2e} ± {sigma_b_rect_fit_1:.1e} m")

# Gráfica vertical
n_fit = np.linspace(min(n_v)-1, max(n_v)+1, 100)

plt.figure()
plt.scatter(n_v, x_v, label="Datos experimentales")
plt.plot(n_fit, A_v_fit*n_fit + B_v_fit, label="Ajuste lineal")

plt.xlabel("m")
plt.ylabel("y (m)")
plt.title("Difracción por abertura rectangular (D = 154.4 cm) dirección vertical")
plt.legend()
plt.grid()
plt.show()



print("\nD2 = 342.8 cm ------------------------------------------------------------")

# Pasamos D a metros D2
D_rect = D_agujero_2 * 1e-2

# Pasamos a metros
dx_h = np.array(d_min_agujero_horizontal_2) * 1e-2
dx_v = np.array(d_min_agujero_vertical_2) * 1e-2

# -----------------------
# Cálculo con Δx D2

A_h = np.mean(dx_h)
A_v = np.mean(dx_v)

sigma_A_h = np.std(dx_h, ddof=1)
sigma_A_v = np.std(dx_v, ddof=1)

# Dimensiones
a_rect_2 = lambd * D_rect / A_h
b_rect_2 = lambd * D_rect / A_v

sigma_a_rect_2 = sqrt(((lambd)/(A_h))**2 * s_D**2 + ((lambd*D_rect)/(A_h**2))**2 * sigma_A_h**2)
sigma_b_rect_2 = sqrt(((lambd)/(A_v))**2 * s_D**2 + ((lambd*D_rect)/(A_v**2))**2 * sigma_A_v**2)

print("\n--- Método Δx ---")
print(f"A_h (para a) = {A_h:.2e} ± {sigma_A_h:.1e} m")
print(f"A_v (para b) = {A_v:.2e} ± {sigma_A_v:.1e} m")
print(f"a = {a_rect_2:.2e} ± {sigma_a_rect_2:.1e} m")
print(f"b = {b_rect_2:.2e} ± {sigma_b_rect_2:.1e} m")

# -----------------------
# RECONSTRUCCIÓN + AJUSTE (horizontal) D2

x_h = np.array([
    - (dx_h[0] + dx_h[1]),
    - dx_h[1],
    + dx_h[2]
])

n_h = np.array(n_agujero_2)

p_h, cov_h = np.polyfit(n_h, x_h, 1, cov=True)
A_h_fit, B_h_fit = p_h
sigma_A_h_fit = np.sqrt(cov_h[0,0])
sigma_B_h_fit = np.sqrt(cov_h[1,1])

a_rect_fit_2 = lambd * D_rect / A_h_fit
sigma_a_rect_fit_2 = a_rect_fit_2 * sigma_A_h_fit / A_h_fit

print("\n--- Ajuste horizontal ---")
print(f"A = {A_h_fit:.2e} ± {sigma_A_h_fit:.1e} m")
print(f"B = {B_h_fit:.2e} ± {sigma_B_h_fit:.1e} m")
print(f"a = {a_rect_fit_2:.2e} ± {sigma_a_rect_fit_2:.1e} m")

# Gráfica horizontal
n_fit = np.linspace(min(n_h)-1, max(n_h)+1, 100)

plt.figure()
plt.scatter(n_h, x_h, label="Datos experimentales")
plt.plot(n_fit, A_h_fit*n_fit + B_h_fit, label="Ajuste lineal")

plt.xlabel("n")
plt.ylabel("x (m)")
plt.title("Difracción por abertura rectangular (D = 342.8 cm) dirección horizontal")
plt.legend()
plt.grid()
plt.show()

# -----------------------
# RECONSTRUCCIÓN + AJUSTE (vertical) D2

x_v = np.array([
    - (dx_v[0] + dx_v[1]),
    - dx_v[1],
    + dx_v[2]
])

n_v = np.array(n_agujero_1)

p_v, cov_v = np.polyfit(n_v, x_v, 1, cov=True)
A_v_fit, B_v_fit = p_v
sigma_A_v_fit = np.sqrt(cov_v[0,0])
sigma_B_v_fit = np.sqrt(cov_v[1,1])

b_rect_fit_2 = lambd * D_rect / A_v_fit
sigma_b_rect_fit_2 = b_rect_fit_2 * sigma_A_v_fit / A_v_fit

print("\n--- Ajuste vertical ---")
print(f"A = {A_v_fit:.2e} ± {sigma_A_v_fit:.1e} m")
print(f"B = {B_v_fit:.2e} ± {sigma_B_v_fit:.1e} m")
print(f"b = {b_rect_fit_2:.2e} ± {sigma_b_rect_fit_2:.1e} m")

# Gráfica vertical
n_fit = np.linspace(min(n_v)-1, max(n_v)+1, 100)

plt.figure()
plt.scatter(n_v, x_v, label="Datos experimentales")
plt.plot(n_fit, A_v_fit*n_fit + B_v_fit, label="Ajuste lineal")

plt.xlabel("m")
plt.ylabel("y (m)")
plt.title("Difracción por abertura rectangular (D = 342.8 cm) dirección vertical")
plt.legend()
plt.grid()
plt.show()


# Presentar las medias de a y de b:

a_media_cuadrada = np.mean([a_rect_1, a_rect_fit_1, a_rect_2, a_rect_fit_2])
s_a_media_cuadrada = np.mean([sigma_a_rect_1, sigma_a_rect_fit_1, sigma_a_rect_2, sigma_a_rect_fit_2])
b_media_cuadrada = np.mean([b_rect_1, b_rect_fit_1, b_rect_2, b_rect_fit_2])
s_b_media_cuadrada = np.mean([sigma_b_rect_1, sigma_b_rect_fit_1, sigma_b_rect_2, sigma_b_rect_fit_2])
print("")
print("-"*75)
print(f"\t\t\t\t\t -----> a = {a_media_cuadrada:.2e} +- {s_a_media_cuadrada:.2e} m")
print(f"\t\t\t\t\t -----> b = {b_media_cuadrada:.2e} +- {s_b_media_cuadrada:.2e} m")

# -----------------------
# Microscopio

rect = np.array(tamano_cuadrado) * 1e-6

a_micro_rect = np.mean(rect[0])
b_micro_rect = np.mean(rect[1])

sigma_a_micro_rect = np.std(rect[0], ddof=1)
sigma_b_micro_rect = np.std(rect[1], ddof=1)

print("\n///Ancho y Alto Agujero por Microscopio///" + "="*33)
print("\n--- Microscopio ---")
print(f"a_micro = {a_micro_rect:.2e} ± {sigma_a_micro_rect:.1e} m")
print(f"b_micro = {b_micro_rect:.2e} ± {sigma_b_micro_rect:.1e} m")




# Porque no coinciden las mediciones del microscopio con las de los ajustes?

# El microscopio nos permite medir el tamaño geométrico, lo que vemos directamente;
# el modelo que usamos se basa en el modelo ondulatorio de la materia, depende de muchas
# aproximaciones. Además la precisión de medición; nosotros no medimos puntos matemáticos,
# medimos zonas oscuras con la imprecisión que eso conlleva de un mínimo no muy bien
# definido. Esto en el orden de mm cuando la medición afecta al orden de micrometros.
# A mayores, nosotros hemos usado la diferencia entre mínimos $\Delta x = x_{n+1}-x_n$
# en lo cual hemos introducido correlaciones entre datos, perdiendo ligera información
# de por medio que amplifica nuestros errores.
# El posible primer fallo es la aproximación del ángulo pequeño de $sin(\theta)=\frac{x}{D}$
# para cuando $D>>x$. Lo cual no deja de ser una aproximación, estamos en el caso de
# un Fraunhofer aproximado con distancias finitas que inducen a deformaciones de imagen.
# A mayores de que el sistema no es ideal, la alineación del láser (si no está centrado
# perfectamente y/o no incide realmente perpendicular al plano) nos afecta al patrón y
# a la medición de la distancia rendija-imagen. Todos estos errores se acumulan y nos trae
# la diferencia teórica-práctica que observamos.